All-electron GW approximation with the mixed basis expansion based on the 

full-potential LMTO method 
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We present a new all-electron, augmented-wave imple- 
mentation of the GW approximation using eigenfunctions 
generated by a recent variant of the full-potential LMTO 
method. The dynamically screened Coulomb interaction W 
is expanded in a mixed basis set which consists of two contri- 
butions, local atom-centered functions confined to muffin-tin 
spheres, and plane waves with the overlap to the local func- 
tions projected out. The former can include any of the core 
states; thus the core and valence states can be treated on an 
equal footing. Systematic studies of semiconductors and insu- 
lators show that the GW fundamental bandgaps consistently 
fall low in comparison to experiment, and also the quasiparti- 
cle levels differ significantly from other, approximate methods, 
in particular those that approximate the core with a pseu- 
dopotential. 
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The GW approximation (GWA) of Hedin 0] has been 
applied to many kinds of materials In customary 

ab initio implementations, the self-energy £ is generated 
from eigenvalues and eigenfunctions calculated within the 
self-consistent local-density approximation (LDA). It has 
been shown that quasiparticle energies computed in this 
way are in significantly better agreement with experiment 
that are the LDA eigenvalues. 

The various implementations of GWA by may classi- 
fied by what kinds of basis sets are used in the expan- 
sion of the LDA eigenfunctions, and the expansion of the 
bare and screened Coulomb interactions v and W. The 
most common implementations make an additional pseu- 
dopotential approximation for the core, which makes it 
possible to expand all these quantities in plane waves. 
However, a plane-wave basis is poorly suited to localized 
orbitals such as d- or /- states. Moreover, as we show 
here, it appears that the pseudopotential approximation 
is somewhat inadequate when used in conjunction with 
the GW approximation. Two other GW implementa- 
tions that do not use pseudopotentials have also been 
published In both of these methods, v and W are 

expanded in plane waves; owing to the difficulty in a 
plane-wave expansion of localized orbitals they did not 
take into account core contributions. Aryasetiawan and 
collaborators implemented a method that expands v and 
W in a linear combination of augmented wave function 
products (product basis) and has applied it to several 
kind of materials, including NiO M, with reasonable re- 



sults. However this implementation requires the atomic 
spheres approximation (ASA) for the LDA, which ap- 
proximates space by a superposition of atom-centered 
("muffin-tin") spheres, neglecting the interstitial. Thus, 
its reliability is uncertain. 

We present a new GW implementation which uses a 
mixed-basis expansion for v and W. v and W are ex- 
panded in Aryasetiawan's product-basis in the muffin-tin 
(MT) spheres, and the interstitial plane waves (IPW) in 
the interstitial region. An IPW is a plane wave with 
the MT contributions projected out. This basis can be 
an efficient one applicable for the localized electrons and 
cores. Together with the rather accurate eigenfunctions 
given by the all-electron full-potential LMTO method, 
our GW method can be both efficient and accurate. It 
is our knowledge the first implementation of GW that 
makes no significant approximations beyond the lack of 
self-consistency (pseudopotential approximation for the 
core, shape approximation in the potential, plasmon-pole 
approximation for the dynamical screening). 

The GW self-energy is 

E(r, r', u>) = ^-J dw'G(r, r',u + uj')e iSuJ 'W(r, r', to'), 

(1) 

where W is the RPA screened Coulomb interaction and G 
is the Green's function. W and G are calculated directly 
from the self-consistent LDA eigenvalues and eigenfunc- 
tions, with exchange-correlation potential V xc {r). Then 
the quasi-particle energy E n (k) (n is the band index and 
k is the wave vector) is given by 

E n (k) - e„(k) = Z nk x 

[(* kn |E(r,r', e „(k))|* kn ) - (* k „|V; L c DA (r)|v|/ k „)] (2) 

where e n (k) and \l/ krl denotes the LDA eigenvalues and 
eigenfunctions, and Z„ k is the quasi-particle (QP) renor- 
malization factor 

Z n]i = [1 - (tf kn |^E(r,r / ,e n (k))|tf kn )]- 1 . (3) 

Our method to calculate the quasi-particle energy, 
is based on the FP-LMTO method which expands the 
eigenfunction ^&k n (we omit the spin index to simplify 
expressions) as 



(4) 
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Here x k (r) denotes the Bloch sum of the MT or- 
bitals Xs( r )- s = ( a i n j^) is a composite index la- 
belling the atom in the unit cell a, the angular mo- 
mentum L, and another index n to specify which MT 
orbitals (multiple envelope functions per L are permit- 
ted) In the interstitial, Xs( r ) ma y be expanded in plane 
waves, while at the augmentation spheres it is matched 
smoothly and onto of a linear combination A au (r) = 
{^aLM^i^), ^ahi^Y^r)} of solutions of the radial 
Schrodinger equation, 4>aL{r), and its energy derivative, 
4>aL(r)- Finally, u = (L,Ip) is a composite of L and ip, 
where Ip assumes value for <\> and 1 for <f>. Therefore 
can be written as 
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where the IPW Pg(i") is defined as 

P k (r) = if r e Any MT 

= exp(i(k + G)r) Otherwise. 



(6) 



In the interstitial, it is evident that products ^kin x v E , k 2 ri' 
can similarly be expanded in IPW's because Pq* (r) x 
Pq 1 (r) = PGi+G 2 ( r )- Within the augmentation sphere 
products of \f are expanded in P k *+ 2 (r), which is the 
Bloch sum of the product basis B am (r). They are con- 
structed from the products of A au (r) x A a / U > (r) following 
the procedure by Aryasetiawan Q. 

Thus a basis suitable for expansion of wave function 
products, and therefore also v and W , is the mixed basis 
{M k (r)} = {Pg(r), P k m (r)}, where Pg is an orthonor- 
mal basis rendered from Pq. I = {G,am} is a compos- 
ite index labelling the basis. The complete information 
to calculate the self-energy and E n (k) are the Coulomb 
matrix u/j(k) = (Mf \v\Mj), the matrix elements of the 
products qn\^ q-im' Mf) and the eigenvalues e^ n . The 
exchange part of the self-energy is written as 
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The screened Coulomb interaction Wu (q, u>) is calcu- 
lated through W = (1 — vD)~ 1 v, where the polarization 
function P/ is 
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77)- (8) 



Finally, the correlation part of self-energy is calculated 
by 



Here the denominator — iS is for occupied states, and +iS 
for unoccupied states. We use the w'-integral method 
given by Aryasetiawan S. Our GW code is developed 
starting from his code. 

Details of the method will be described elsewhere. Ap- 
plying the method to a wide range of semiconductors, we 
find a number of systematic tendencies. The most impor- 
tant ones are: (i) there is a systematic underestimate of 
the fundamental bandgap (it) The bandgap error tends 
to increase with bond polarity, and for compounds with 
relatively shallow cation d states. (Hi) the GW rather 
badly underestimates the deepening (relative to LDA) of 
occupied cation d states found in III-V and II- V com- 
pounds, (iv) For zincblende lattices, the bandgap error 
tends to be slightly larger at the X point than at V. (v) 
the bandgap is sensitive to what approximations are used 
for shallow core states. 

GaAs illustrates many of these general findings. It has 
a direct gap at OK of 1.52 eV; however to compare against 
the present non spin-polarized calculations, we assign the 
valence band maximum r 15 to + |r 8u . (Spin-orbit 
coupling splits r 15 by shifting the r 8 „ states +^A and 
the Tr v state — |Ao. This was confirmed by an LDA- 
ASA calculation including spin-orbit coupling, which also 
yielded Ao=0.36 eV, in good agreement with the ob- 
served Ao=0.34 eV.) Table | shows gaps at T, X and 
L; the spin-orbit-corrected fundamental gap is 1.63 eV. 
An accurate LDA calculation must include both 3d states 
(which push upwards on Tis, narrowing the gap) and the 
Ad states (which push downwards, widening the gap). At 
present the GW implementation can include only one of 
these states for the input wave functions (the GW it- 
self includes both) , but to compute the bands within our 
present GW we must select either 3d or Ad to input to 
the GW. Therefore, all three LDA cases are presented 
(the 3+ Ad case uses local orbitals pQ]), and are seen to 
fall within ^0.15 eV of each other. 

Similarly the GW gaps, whether the LDA 3c? basis or 
Ad basis is input, lie within ~0.15 eV of each other (albeit 
with a slightly different /c-dependence) . By comparing 
LDA bands to the 3+Ad case, the 4c? is a better choice for 
GaAs, especially in the conduction bands where the Ga 
4c? begin to play an important role. In any case, the di- 
rect gap (rx5„ — > Tic) is ^0.3 eV smaller than the exper- 
imental value; while the X and L point (Xi c and Xi c ) are 
underestimated by a somewhat larger amount, ~0.4 eV. 
We have found this tendency to be rather systematically 
followed in the III-V semiconductors, as shown below. 

The effect of the core is particularly important. Line 
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coreO in Table Q shows what gaps result when valence 
states only are included in the calculation of D and £. 
For this case the LDA potential we subtract corresponds 
to the valence-only density, V^P (n va i ; r) . This approx- 
imates what is typically done in pseudopotential GW 
(GW:PP) calculations. It is seen that this approxima- 
tion leads to much too strong a A-dependence on the gap 
shift, something also seen when compared to GW:PP cal- 
culations. Table | shows some data of Shirley et. al. [ Q |, 
(marked GW:PP). They rather closely track the coreO 
results except for an approximately fc-independent shift 
of 0.32 eV. Shirley also included an approximate core 
polarization term (marked GW:PP+CP). It is seen that 
the addition of core-polarization terms to the PP have 
roughly similar effects (increasing the shift at T, decreas- 
ing it at X), but the shifts are larger than the all-electron 
results, and depend more strongly on k. Similar overes- 
timates of gaps by the GW:PP method are seen in Si 
and AlAs. For example Shirley computed the Ti$ v ^Xi c 
transition to be 1.31 eV, while our all-electron result is 
1.04 eV, close to the GW:PAW result (1.10 eV) of Arnaud 
[ 5|, and an early GW:LAPW result (1.14 eV) by Hamada 
[1|. Moreover, we find for GaN (and indeed generally for 
all semiconductors studied), a weak energy-dependence 
of the GW-LDA shift of the conduction bands up to 10 
eV above the conduction-band minimum, in accord with 
a careful analysis of both UV reflectivity and near-edge x- 
ray absorption spectra ll2] . However, using the GW:PP 
method, Rubio et. al ]13| ] found a significant energy- 
dependence of this shift. 

Inspecting the bottom of Table GL it is evident that 
the Ga 3d states contribute in an important way to the 
dielectric response. The row marked core2 shows that 
the levels change rather significantly when the Ga 3d is 
omitted from the calculation of D. 

Fig. [I] illustrates the systematics of the errors in the 
fundamental gap for a range of semiconductors and in- 
sulators. With the sole exception of diamond [14], the 
GW gaps are systematically smaller than experiment. 
As noted for GaAs above, there is some fc-dependencc 
to the errors: the lowest conduction band levels at L and 
more especially at X are generally in worse agreement 
than those at T. Two sources of error are readily identi- 
fiable, and as will be shown elsewhere |lq| , they account 
for nearly all of the gap errors in these sp bonded semi- 
conductors. As noted by Maksimov []l6|| , an important 
contribution to the gap error arises from the nonlocal 
screened exchange missing from LDA. In Hartree-Fock 
theory, the bare exchange is used, therefore it rather 
severely overestimates the bandgap and bandwidths in 
general. The GW does screen the exchange, and qualita- 
tively speaking GW replaces in addition to V^ DA , a term 
which corresponds to the nonlocal screened exchange, 
with the nonlocal screened exchange written as 



K L c DA + -(Sx 



V, 



LDAn 



(10) 



In the present non self-consistent GW implementation, 
e is computed from the LDA, with e w ^10 to 20% over- 
estimated on account of the LDA gaps being too small. 
As we will show elsewhere |fL5fl , self-consistency reduces 
e slightly. This enhances the nonlocal part of the self- 
energy operator, and further increases the GW-LDA gap 
correction, as can be qualtitatively understood by Eq. [l^. 
This correction is 0.2 eV in Si |l7j but increases in the 
more ionic materials where e is smaller. That the cor- 
rection is larger in the latter case might be expected 
because the initial LDA potential (from which GW is 
constructed) is a poorer approximation, as again can be 
qualitatively seen from Eq. |l0|. 

Table || shows the underbinding by GW of the oc- 
cupied cation d levels. As is shown, the LDA rather 
badly underestimates this binding. The LDA error is well 
known and also reasonably well understood Jl8[ ] , namely 
that the LDA eigenvalues are inappropriately interpreted 
as excitation energies. It is significant that, while the 
GW shifts these levels in the direction of experiment, the 
shift is much too small. Comparison with the exchange- 
only (Hartree-Fock) column shows clearly that the GW is 
strongly overscreening the bare exchange for these cation 
d levels. (Similar difficulties are found for the so-called 
charge transfer antiferromagnets such as NiO and MnO, 
as will be discussed elsewhere.) This underbinding also 
contributes an important term to the gap error; shows 
importance increases for shallow d states. It is not clear 
whether self-consistency in the GW will remedy this er- 
ror. It is interesting that, even within the LDA, these lev- 
els can be reasonably computed using Slater transition- 
state theory. 

The underbinding of the d levels makes an important 
contribution to the gap error where the levels are rela- 
tively shallow (e.g. CdTe, and GaAs). There is a cou- 
pling between the d state and and the valence band max- 
imum, which pushes the latter upwards and reduces the 
bandgap. When the energy separation between these 
states are too small, the coupling them is overestimated. 
This will be shown in some detail in a future work [ fl5| |. 
In additon, we will publish the GW results on the wurzite 
ZnO where we prepare the input cigenfunctions by the 
full-potential LAPW method [||. 

In summary, an all-electron implementation of GW 
presented here leads to systematic errors in semiconduc- 
tor and insulator bandgaps, whose origin were identified. 
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TABLE I. Selected energy eigenvalues, in eV, at F, L and 
X for GaAs. Experimental data are taken from Ref. Jl9| ; the 
line "Expt-SO" subtracts the spin-orbit coupling from the raw 
experimental data to compare with the calculated levels, as 
discussed in the text. The GW:PP+CP and GW:PP data 
are GW results jll] using an LDA pseudopotential, with and 
without a core polarization term added. The GW results are 
quoted with the Z (underlined), and with Z—l immediately 
below. Also shown are changes in the GW levels owing to 
varying degrees of approximation for the Ga and As 3d states. 





Tic 


Tl5c 


Lie 


Xl c 


X3 C 


Expt 


1.52 


4.72 


1.84 


2.01 


2.41 


Expt-SO 


1.63 


4.83 


1.95 


2.12 


2.52 


t n a 


u.oo 


o.uo 


u.ou 


1 "id 




LDA(4d) 


0.41 


3.72 


0.91 


1.34 


1.58 


LDA(3d) 


0.24 


3.57 


0.76 


1.33 


1.57 


GW,4d(Z) 


1.35 


4.20 


1.54 


1.56 


1.97 


GW,4d(Z=l) 


1.50 


4.28 


1.63 


1.59 


2.02 


GW,3d(Z) 
GW,coreO d 


1.34 


4.38 


1.59 


1.71 


2.04 


0.94 


4.18 


1.37 


1.72 


2.01 


GW:PP 


1.29 




1.69 


2.05 


2.34 


GW:PP+CP 


1.69 




1.88 


1.93 


2.29 


HF,4d 


4.84 


8.50 


5.10 


5.26 


5.63 


<5GW,core3 a 


0.02 


0.03 


0.02 


0.03 


0.04 


<5GW,core2 b 


0.24 


0.34 


0.26 


0.25 


0.24 


<5GW,corel c 


0.26 


0.35 


0.27 


0.25 


0.25 


5GW,coreO d 


-0.41 


-0.02 


-0.17 


0.16 


0.04 



a core3 treats the As 3d at the Hartree-Fock level (bare ex- 
change) and omits this state in the computation of D 
h core2 additionally omits the Ga 3d in the calculation of D 
c corel additionally approximates the Ga 3d at the Hartree- 
Fock level 

d coreO computes both D and W with valence states only, com- 
pletely neglecting the core beyond the LDA treatment. 



TABLE II. Occupied cation d band energy levels in se- 
lected II- VI and III-V semiconductors, relative to the va- 
lence band maximum. Column 'HF' is the exchange-only 
(Hartree-Fock) result. 





Expt 


LDA 


GW 


HF 


GaN 


-17.1 a 


-13.6 


-16.4 


-23.3 


GaAs 


-18.8 C 


-15.0 


-18.1 


-27.2 


InP 


-16.8 C 


-14.2 


-15.7 


-21.2 


InAs 


-17.1 c 


-14.4 


-16.1 


-22.6 


ZnS 


-8.7 b 


-6.2 


-7.1 


-14.1 


ZnSe 


-9.0 C 


-6.7 


-7.7 


-15.2 


CdS(ZB) 


-9.2 d 


-7.5 


-8.2 


-12.3 


CdTe 


-10.5 C 


-8.1 


-9.3 


-14.4 


HgTe 


-8.6 C 


-7.1 


-7.6 


-12.8 
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average of X and F, "turning points" method. 
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Error in calculated bandgap 




- GW error, direct gap 
GW error, fundamental gap 
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C AlSb GaSb InSb CdTe ZnS MgO 
Si AlAs GaAs InAs HgTe GaN 
Ge AIP GaP InP ^nSe InN 

FIG. 1. Top panel: experimental values of the fundamental 

gap in a variety of semiconductors and insulators, approxi- 
mately ordered by ionic character. Bottom panel: deviations 
from experiment in the fundamental gap. Red points are LDA 
errors; green points are GW errors. The light and green de- 
notes the smallest direct gap; the dark green denotes the fun- 
damental gap when it differs from the direct gap. A term 
1/3 Ao was subtracted from the calculated GW and LDA data 
to account for spin-orbit coupling, as discussed in the text. 
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